// LPMolecule.cpp: implementation of the LPMolecule class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "lpmolecviewercontrol.h"
#include "LPMolecule.h"

//#include "StringTokenizer.h"
#include "math.h"
#include <new>

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
//#define new DEBUG_NEW
#endif

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

LPConnector *lcon;

LPMolecule::LPMolecule()
{
//Constructor
	this->CONNECTOR_MAX_NUMBER=5000;
	this->nconnectors=0;
	this->ATOM_MAX_NUMBER=3000;
	this->natoms=0;
	ListAtoms=new LPAtom*[this->ATOM_MAX_NUMBER];
	ListConnectors=new LPConnector*[this->CONNECTOR_MAX_NUMBER];
	AUTOCONNECTDISTANCE=1.85; // Default autoconnect distance in angstrons
	//ctoken=NULL;
}



LPMolecule::~LPMolecule()
{
	for (int i=0; i<natoms;i++){
		delete ListAtoms[i];
	}
	delete[] ListAtoms;
	
	for (i=0; i<nconnectors;i++){
		delete ListConnectors[i];
	}
	delete[] ListConnectors;

}

int LPMolecule::LoadXYZMolecule(char *fp)
{
	
	//Reads file pointer fp and sets up molecule structure accordingly
	char Line[200];
	char Name[12];
    double xpos, ypos, zpos;
    int atoms;
	double charge, u, v, w; //Not used yet
 
    //register Atom __far *ptr;  //Atom is a struct defined in molecule.h
	LPAtom *ptr1;
    int i,count;
 
	//Tokenizes the file
	//StringTokenizer *tokenizer=new StringTokenizer(DataFile,'\n');
	char *tokenizer=strtok(fp,"\n");

	// Read Number of Atoms
    //strcpy(Record,tokenizer->getNextElement());
    strcpy(Line,tokenizer);

	sscanf(Line,"%d",&atoms); //Reads number of atoms and stores in &atoms
	if (atoms>1 ||atoms <1000){

		tokenizer=strtok(NULL,"\n"); //Read next line with molecule name
		if (tokenizer!=NULL) (strcpy(Line,tokenizer));

		strcpy(this->InfoMoleculeName,Line);
		//Reads molecule name and place it in InfoMoleculeName variable
		for( i=0; i<atoms && tokenizer!=NULL; i++ )
        {   
		//strcpy(Record,tokenizer->getNextElement());
            tokenizer=strtok(NULL,"\n");
			strcpy(Line,tokenizer);
			ptr1=new LPAtom();
            ptr1->serno=i;
			
			//Line must contain  data
            xpos = ypos = zpos = 0.0;
            count = sscanf(Line,"%s %lf %lf %lf %lf %lf %lf %lf",
                           Name, &xpos, &ypos, &zpos, &charge, &u, &v, &w );
			//Reads description of atom
 
 
			//TODO:From atom type process elem number
			ptr1->setType(Name);
			//strcpy(ptr1->type,type);

			//Stores position description of atom
			ptr1->x = xpos;
            ptr1->y = ypos;
            ptr1->z = zpos;      

			this->ProcessAtom(ptr1,m_autoconnect); //Adds atom to molecule
			//ProcessAtom( ptr );
        }
    }

	return( true );
}



void LPMolecule::ProcessAtom(LPAtom *ptr,bool autoconnect)
{
	double dist;
	//Adds atom to list
	if (this->natoms<ATOM_MAX_NUMBER){
		ListAtoms[this->natoms]=ptr;
		natoms++;
	}

	//TODO: Create Connectors? automatically.
	if (autoconnect==true){
		//Calculates distance from this atom to the others
		//If distance is bellow a certain value create conector
		LPAtom *ai;
		for (int i=0;i<this->natoms-1;i++){
			ai=ListAtoms[i];
			dist=sqrt((ai->x-ptr->x)*(ai->x-ptr->x)+
					(ai->y-ptr->y)*(ai->y-ptr->y)+
					(ai->z-ptr->z)*(ai->z-ptr->z));
			if (dist<AUTOCONNECTDISTANCE){
			//Forms connection	
				//lcon=new LPConnector();
				//lcon->Atom1=ai;
				//lcon->Atom2=ptr;

				//Connects here, but
				//does not connect if both elements are hydrogen atoms
				if (ai->elemno!=1 || ptr->elemno!=1) {  
					ListConnectors[nconnectors]=new LPConnector();
					ListConnectors[nconnectors]->Atom1=ai;
					ListConnectors[nconnectors]->Atom2=ptr;
					nconnectors++;
				}
				
			}
		}

	}
};

int LPMolecule::getNumberAtoms()
{
	return(this->natoms);
}

void LPMolecule::centerMolecule()
{
if (this->natoms>0){
	int i;
	double xsum=0,ysum=0,zsum=0;
	double xav,yav,zav;
	
	for (i=0;i<this->natoms;i++){
		xsum+=this->ListAtoms[i]->x;
		ysum+=this->ListAtoms[i]->y;
		zsum+=this->ListAtoms[i]->z;
	}
	xav=xsum/natoms;
	yav=ysum/natoms;
	zav=zsum/natoms;	
	for (i=0;i<this->natoms;i++){
		xsum=this->ListAtoms[i]->x-xav;
		ysum=this->ListAtoms[i]->y-yav;
		zsum=this->ListAtoms[i]->z-zav;

		this->ListAtoms[i]->x=xsum;
		this->ListAtoms[i]->y=ysum;
		this->ListAtoms[i]->z=zsum;
	}


}	
}

double LPMolecule::getDimension()
{
	double max=0;
	if (natoms>0){
		for (int i=0;i<natoms;i++){
			if (max<fabs(ListAtoms[i]->x)) max=fabs(ListAtoms[i]->x);
			if (max<fabs(ListAtoms[i]->y)) max=fabs(ListAtoms[i]->y);
			if (max<fabs(ListAtoms[i]->z)) max=fabs(ListAtoms[i]->z);
		}
	}
	return max;
}

int LPMolecule::getNumberConnectors()
{

	return(this->nconnectors);

}

bool LPMolecule::LoadMolecule(TCHAR *filename, bool autoconnect)
{
	long filesize;
	char *filechar;
	
	m_autoconnect=autoconnect;

	//FILE *file=fopen(filename,"rb");
	FILE *file=_tfopen(filename,"rb"); //TCHAR (unicode) version
	
	if (file!=NULL){
		int i=fseek(file,0,SEEK_END);
		filesize=ftell(file);
		filechar=new char[filesize+1];
		rewind(file);
		i=fread(filechar,1,filesize,file); //reads the whole file
		filechar[filesize]='\0';	//null terminates string
		fclose(file);				//Close file
		string ext=getFileExt(filename);
		if (strcmp(ext.c_str(),"xyz")==0 || strcmp(ext.c_str(),"XYZ")==0){
			this->LoadXYZMolecule(filechar);
			return (true);
		}
		if (strcmp(ext.c_str(),"pdb")==0 || strcmp(ext.c_str(),"PDB")==0){
			this->LoadPDBMolecule(filechar);
			return (true);
		}

	}
	return (false);
}

string LPMolecule::getFileExt(const string &s)
{
   size_t i = s.rfind('.', s.length( ));
   if (i != string::npos) {
      return(s.substr(i+1, s.length( ) - i));
   }

   return("");

}

int LPMolecule::LoadPDBMolecule(char *fp)
{
	//Reads file pointer fp and sets up molecule structure accordingly
	char ResName[3]="";
	char Record[7]="";
	char value[8]="";
	char Line[200];
	char *ctemp1,*ctemp2;

    double xpos, ypos, zpos;
    //int atoms; //number of atoms
	//double charge, u, v, w; //Not used yet

	LPAtom *ptr1;
    int i;
 
	i=0;
	//Tokenizes the file for reading lines
	char *tokenizer=strtok(fp,"\r\n");

	while (tokenizer!=NULL){
		//Line-by-line processing

		int ilen=strlen(tokenizer);
		strncpy(Line,tokenizer,sizeof(Line)); //Copies to variable Line
		
		//Check of occurrence of END or TER
		ctemp1=strstr(Line,"END");
		ctemp2=strstr(Line,"TER");
		if (ctemp1==NULL && ctemp2==NULL){
			//String not found, Can continue...
		
			if (ilen>6){
				ctemp1=strstr(Line,"HETATM");
				ctemp2=strstr(Line,"ATOM");
				if ((ctemp1!=NULL || ctemp2!=NULL) && ilen>53){
					//Atom found
					ptr1=new LPAtom();
					ptr1->serno=i;
					strncpy(ResName,Line+12,3); //Copies 3 characters, note strncpy does not null-terminate
					ResName[3]='\0';			//null terminate string
					strncpy(value,Line+31,8);
					value[7]='\0';
					xpos=atof(value);
					strncpy(value,Line+39,8);
					value[7]='\0';
					ypos=atof(value);
					strncpy(value,Line+47,8);
					value[7]='\0';
					zpos=atof(value);

					ptr1->setType(ResName);
					ptr1->x = xpos;
					ptr1->y = ypos;
					ptr1->z = zpos;      

					this->ProcessAtom( ptr1 , false); //Adds atom to molecule
												// but don't add connections yet
					i++;
				}
			
				ctemp1=strstr(Line,"COMPND");
				if (ctemp1!=NULL){
					strncpy(this->InfoMoleculeName,Line+6,sizeof(InfoMoleculeName)); //Copies the name of the atom
					//this->InfoMoleculeName[50]='\0';		//Null-terminate
				}

				ctemp1=strstr(Line,"CONECT");
				if (ctemp1!=NULL){
					LPAtom *ai;
					int n0, n1;
					strncpy(value,Line+7,5); //Reads serial number of central atom
					value[4]='\0';
					n0=atoi(value);
					int slen=strlen(Line);
					for(int i=0; i<4 && slen> 15+(i*5); i++){
						strncpy(value,Line+12+(i*5),5); //Reads serial number of bonded atom
						value[4]='\0';
						n1=atoi(value);
						if (n1>n0 && n0>0){ //Adds connector, avoids repeating connectors
							ListConnectors[nconnectors]=new LPConnector();
							ai=ListAtoms[n0-1];
							ListConnectors[nconnectors]->Atom1=ai;
							ai=ListAtoms[n1-1];
							ListConnectors[nconnectors]->Atom2=ai;
							nconnectors++;
							//Attention, this routine assumes that the serial number
							//of atoms appear oredered when the file is read
						}
					} //for
				} //if
			}
			tokenizer=strtok(NULL,"\r\n");		//Gets new line
		}
		else{	//END or TER found
				//exit while routine by flagging with tokenizer
			tokenizer=NULL;
		}
	}

	//If no connector information in file
	//then add connections here if autonnect=true
	if (m_autoconnect && nconnectors==0){
		connectMoleculeAtoms();
	}
	//delete Line;
	return (true);
	
}

void LPMolecule::connectMoleculeAtoms()
{
	double dist;
	LPAtom *ai;
	LPAtom *aj;
	for (int i=0;i<this->natoms-1;i++){
		ai=ListAtoms[i];
		for (int j=i+1; j<this->natoms;j++){
			aj=ListAtoms[j];
			
			//Does not connect if both elements are hydrogen atoms
			if (ai->elemno!=1 || aj->elemno!=1){

				//Calculate distance between atoms
				dist=sqrt((ai->x-aj->x)*(ai->x-aj->x)+
					(ai->y-aj->y)*(ai->y-aj->y)+
					(ai->z-aj->z)*(ai->z-aj->z));

				if (dist<AUTOCONNECTDISTANCE){
				//Forms connection	
					ListConnectors[nconnectors]=new LPConnector();
					ListConnectors[nconnectors]->Atom1=ai;
					ListConnectors[nconnectors]->Atom2=aj;
					nconnectors++;
				}
			}	
		}
	}
}
